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We report a Monte Carlo finite-size scaling study of the demixing transition of a symmetrical 
Lennard- Jones binary fluid mixture. For equal concentration of species, and for a choice of the 
unlike-to-like interaction ratio 8 = 0.7, this transition is found to be continuous at liquid-vapor 
coexistence. The associated critical end point exhibits Ising-like universality. These findings confirm 
those of earlier smaller scale simulation studies of the same model, but contradict the findings of 
recent integral equation and hierarchical reference theory investigations. 
PACS numbers: 64.70Dv, 02.70.Lq 



I. INTRODUCTION AND BACKGROUND 



Symmetrical binary fluid mixtures are two-component fluid models in which the configurational energy is invariant 
with respect to interchange of the two particle species. One example is the symmetrical Lennard- Jones (LJ) mixture 
in which the interactions between particles of species 1 and 2 are controlled by a LJ potential with scale parameters 
(Tn = C22 = (Ti2 = o~ and interaction strengths en = 622 = e ^ £12 . The phase diagram of the mixture is uniquely 
specified by the ratio of interaction strengths between the unlike and like species: 5 = £12/ e. 

Symmetrical mixtures have been the subject of considerable recent attention on account of their surprisingly rich 
phase behaviour [JHHSOI. The phase diagram of the mixture is typically represented in terms of the temperature 
T, the overall number density p = p± + P2, and the concentration c = p2/(pi + P2)- Most attention has focused on 
the case of equal species concentration, c = 0.5. Within this particular symmetry plane of the full phase diagram, 
the p-T dependence of the phase behaviour on the value of 6 has been systematically studied using simulation and 
mean field theory for the case of a symmetrical square- well mixture jf| . The findings of that work are summarized 
in fig 1. of ref. |5j. Three distinct classes of phase diagram were found depending on the choice of the parameter 5. 
Specifically, for large 6 < 1 (fig. 1(a) of ref. |5j) there occurs a "A-line" of critical demixing transitions which intersects 
(and is truncated by) the liquid- vapor line at a critical end point (CEP). For small 5 > 0, on the other hand, the line 
of demixing transitions intersects the liquid- vapor line at the liquid- vapor critical critical point, forming a tricritical 
point (fig. 1(c) of ref. @). Intermediate between these two regimes, one observes both a liquid-vapor critical point 
and a tricritical point at higher densities (fig. 1(b) of ref. 5]). In the simulations this latter regime was observed to 
occur for 0.65 < S < 0.68, while within the particular mean field theory employed it occurred for 0.605 < <5 < 0.708. 

It has been suggested |5j that the scenario of phase behaviour outlined above for the square well mixture, might 
in fact be generic to all symmetrical fluids. On the theoretical side, support for this proposal has come from studies 
of variety of model systems, including mean field studies of a lattice fluid model @ and integral equation theories 
La m U \1D\ of the hard core Yukawa fluid mixture (HCYFM). The same overall scenario was initially reported for 
the hard core + Lennard- Jones mixture by Pini et al 1 1 1 1 ] in a study employing the Hierarchical Reference Theory 
(HRT)- a powerful fluctuation-inclusive approach based on renormalization group techniques [hi . By contrast, the 
available simulation data for models other than the square well fluid is less comprehensive and has,in the main part, 
concentrated on the regime of large 5. Specifically, MC simulations of the symmetrical LJ fluid |J| found evidence 
for a CEP (and associated coexistence curve singularities) for the particular case S — 0.7. Similarly a CEP has been 
reported in simulation studies of the HCYFM for 6 = 0.9 [Tof . 

Very recently, however, the situation with specific regard to the CEP regime has become less clear following two 
theoretical studies which question its existence in the expected range of S. In the first of these, Antonevych et al. 01 
applied the modified hypcrncttcd chain (MHNC) integral equation theory to the symmetrical LJ fluid, but found no 
evidence for a CEP at c = 0.5. Instead they reported that the intermediate regime (fig. 1(b) of ref. [j|) persists right 
up to 8 = 0.81, (the limit of their study). This value of 5 is considerably greater than that (d = 0.7) at which CEP 
behaviour was originally found in the simulations of the same model. On this basis, Antonevych et al. suggested 
that the demixing transition could always be first order at liquid-vapor coexistence, although they speculated that it 
might be only weakly so. Furthermore, they questioned the ability of simulation to distinguish this possibility from 
a true CEP. 

In a separate study, Pini et al. 0] have applied the HRT approach to the HCYFM. They too observed no CEP, 
instead reporting (in accord with the finding of Antonevych et al. |13|) that the intermediate regime persists (albeit 
weakly) to <J = 0.8, the limit of their study. Moreover, these authors report that in their previous study of the 
LJ mixture the apparent observation of a CEP was an artifact stemming from the low resolution of the grid in 
the p-c plane on which the HRT equations were solved. 
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Clearly, therefore, the evidence emerging from the various theoretical and simulation studies is contradictory with 
regard to the existence of a CEP and the range of values of 8 for which it occurs. Most significant in our view is the 
failure of HRT to observe a CEP for S < 0.8. Given the fluctuation-inclusive nature of this approach (and its success 
in other contexts, se e e.g . |l6|), one cannot simply ascribe discrepancies between the simulations and the MHNC 
study of the LJ fluid [l3| to the neglect of critical fluctuations by the latter. 

Questions have been raised too regarding the reliability of simulation in distinguishing a CEP from a weak first 
order transition |l3j . While, in our opinion, the original work £j does provide compelling evidence for CEP behaviour, 
it is certainly the case that owing to the prevailing computational constraints, the range of system sizes studied was 
smaller than one might have wished. This in turn precluded an investigation of corrections to scaling and the approach 
to the limiting critical behaviour. 

In the light of the above considerations, it would seem worthwhile to revisit the original system in which a CEP was 
observed, with a view to performing a more comprehensive determination of the nature of the demixing transition 
at liquid-vapor coexistence. To this end we have carried out a detailed Monte Carlo (MC) finite-size scaling (FSS) 
simulation study of the symmetrical LJ fluid with 6 = 0.7, for which it was previously claimed that a CEP occurs 
0. The new study has been executed along similar lines to the original one, but employed a considerably larger 
range of systems sizes, thus permitting a more sophisticated finite-size scaling analysis. The results demonstrate 
(unambiguously, we believe) that the model does indeed exhibit a CEP of the expected Ising type. 



II. METHOD AND RESULTS 



The simulation methodology employed here is broadly similar to that described in ref. |J] and we refer the reader 
to that paper for a detailed description. Briefly, we have performed a MC simulation study of the liquid-vapor 
coexistence line of the symmetrical LJ fluid mixture with S = 0.7 in the neighborhood of the demixing transition. 
The simulations were carried out within the grand canonical ensemble and the chemical potentials of the two species 
were constrained to be equal (/ii = fJ-2 = p), implying on symmetry grounds that (c) = 0.5. The system was confined 
to a cubic box of linear dimension L having periodic boundary conditions. In all, seven system sizes were studied 
having L — 10c, 12.5<r, 15<r, 17.5a, 20a, 25a and 30er. The largest of these systems contains some 1.6 x 10 4 particles 
in the liquid phase at the demixing point. The range of system sizes studied here is to be compared with that of 
the original study for which the largest size attained was L — 12.5a. The greater range of sizes studied here was 
partly facilitated by the use of high performance parallel computers. 

The liquid- vapor coexistence curve /x cx (T) was obtained using multicanonical extended sampling techniques |20| . 
These allow both the liquid and vapor phases to be sampled (and thus connected) in a single simulation run. A 
suitable form for the requisite preweighting function was obtained using a variant of the recently proposed Wang- 
Landau method |17| . To locate the coexistence line , th e equal peak area criterion was applied to the measured form of 
the number density distribution function pl{p) |l8U19j . Histogram reweighting techniques were employed to fine-tune 
this procedure and to facilitate the mapping of a portion of the liquid-vapor coexistence curve in the neighborhood 
of the demixing transition. 

Fig. n shows the measured form of the liquid- vapor coexistence line /z cx (T) |2l|. To elucidate the nature of the 
demixing transition along the liquid branch of this tightly determined coexistence line, we have measured the prob- 
ability distribution function p^(m) of the demixing order parameter m, where m = (Ni — A^)/(^Vi + N2), with N\ 
and N2 the instantaneous counts of the respective particle species. The form of pi(m), for the L — 25a system at 
temperatures spanning the demixing transition, is shown in fig. El W e nn d that pi(m) evolves smoothly between 
the strongly double peaked form (corresponding to a demixed liquid) at low temperature and a singly peaked form 
(corresponding to a mixed liquid) at higher temperatures. 

To quantify this evolution more precisely, we have measured the temperature dependence of the fourth order 
cumulant ratio: Ul = 1 — (m 4 ) /3(to 2 ) 2 , which provides a dimensionless measure of the shape of a distribution. The 
results (fig demonstrate (for all system sizes) a smooth evolution of the concentration distribution from the ordered 
to the disordered phase. The rate of change increases with the system size, but there is no evidence of a jump 
discontinuity. More significantly, there is a fairly well defined crossing point of the curves for Ul ~ 0.46 close to the 
reduced temperature T* = 1/c = 0.959. Such an intersection reflects scale invariance in the concentration distribution 
and does not occur at a first order phase transition. We note further that the magnitude of the cumulant ratio at the 
intersection point corresponds closely to the known Ising univer sal value, U L sing = 0.465 [23. 

Although the results of fig.EJconstitute strong evidence for a critical demixing transition, a close inspection reveals 
that the cumulant crossings do not occur precisely at a unique temperature for all system sizes. Such behaviour 
mirrors that observed in other contexts (see eg. [r| 0|) where it was found to be attributable to corrections to 
finite-size scaling. Since these corrections have a known universal scaling form, one can attempt to fit them, thereby 
permitting an extrapolation to the thermodynamic limit. The procedure involves determining, for each system size, 
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the coexistence state point for which the form of pz,(m) on the liquid branch best matches the (independently known 
\22\) fixed point order parameter distribution function p* (to) appropriate to the Ising universality class. The matching 
temperature as a function of system size is expected to scale like T C {L) ~ L~( 1 + e )/ 1/ ^<|, where 9 w 0.50(2) is the 
Ising correction-to-scaling exponent [2j| and v = 0.6294(4) is the correlation length exponent [24|. Fig. 0] shows the 
data expressed in this form, which indeed exhibits the anticipated scaling behaviour. Extrapolation of a linear fit to 
L = oo yields an estimate for the critical temperature T ccp = 0.9595(3), which occurs for ^ cop = —3.704(2). 

It is instructive to plot the finite-size forms of pi, (to) at the estimated critical point in order to expose the scale 
and character of corrections to finite-size scaling. Figure compares the distributions for the L = 12.5ct and L = 25a 
system sizes with that of p*(m) |22J. One sees that the results for the larger system size closely match the limiting 
form. 

Finally, we have considered the behaviour of the number density along the liquid branch at coexistence. Were the 
demixing transition first order along the coexistence line, one would expect that coupling between the concentration 
and the number density would engender a first order phase transition in the liquid branch number density (cf. fig. 
1(b) of ref. f5j). This would be manifest as a discontinuity in the density of the coexisting liquid phase and by a 
splitting of the liquid-phase peak in the density distribution pl{p). Fig. shows the measured form of the liquid-phase 
peak for the L = 25<r system size, for several coexistence state points spanning the demixing temperature. Also shown 
(inset) is the average liquid phase density. The distribution of the density in the liquid phase is singly peaked and we 
find no evidence of a discontinuity in the average density for any of the system sizes studied |25| . 



III. DISCUSSION AND CONCLUSIONS 



In summary, we have performed a finite-size scaling analysis of the demixing transition of a symmetrical LJ fluid 
along the liquid branch of the liquid-vapor coexistence line. The results demonstrate that for the case 5 = 0.7, this 
transition is continuous and of the expected Ising type 26]. Accordingly we can conclude that this system does indeed 
exhibit a critical end point, as originally proposed in rcf 4J. This clear finding should provide a benchmark for testing 
the reliability of theoretical approaches to the phase behaviour of symmetrical mixtures. 

With regard to the original motivation for this study, namely the failure of both the MHNC integral equation study 
and HRT to observe a CEP in the expected range of 6 13, Il4| . it is difficult to comment definitively at this stage 
on the source of the discrepancies. We note, however, that problems with the accuracy of the MHNC approach for 
studying the liquid branch of the HCYFM have previously been reported . With regard to the HRT, it would 
certainly be useful if this method could be applied to other symmetrical fluid models in order to determine whether 
there is any systematic model dependence on the qualitative features of the phase behaviour. 
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FIG. 1: The liquid-vapor coexistence curve in the neighborhood of the demixing transition, obtained for the L = 25<r system. 
Statistical errors are considerably smaller than the symbol sizes. 




FIG. 2: The distribution of the order parameter obtained for the L = 25a system at a number of temperatures along the liquid 
branch of the coexistence curve. In each case, the order parameter has been scaled m-tmto ensure that all distributions have 
unit variance. 
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FIG. 3: The temperature dependence of the fourth order cumulant ratio (see text) measured along the liquid branch of the 
coexistence curve. 




FIG. 4: The finite-size scaling behaviour of the corrections to scaling obtained according to the method described in the text 
and ref. fl9|l . The straight line is a linear fit to the data. 
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FIG. 5: Order parameter distribution for the L = 12.5a and L = 25a system sizes at the extrapolated parameters of the critical 
end point (cf. figHJ. Also shown (solid line) is the limiting universal fixed point form p*(rh) |2^|. 




FIG. 6: The measured form of Pl(p) on the liquid branch for the L = 25a system size for a number of temperature spanning 
the demixing point, The inset shows the corresponding average number density as a function of temperature. 



